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Abstract 

Background: In order to reduce time and efforts to develop microbial strains with better capability of producing 
desired bioproducts, genome-scale metabolic simulations have proven useful in identifying gene knockout and 
amplification targets. Constraints-based flux analysis has successfully been employed for such simulation, but is 
limited in its ability to properly describe the complex nature of biological systems. Gene knockout simulations are 
relatively straightforward to implement, simply by constraining the flux values of the target reaction to zero, but the 
identification of reliable gene amplification targets is rather difficult. Here, we report a new algorithm which 
incorporates physiological data into a model to improve the model's prediction capabilities and to capitalize on the 
relationships between genes and metabolic fluxes. 

Results: We developed an algorithm, flux variability scanning based on enforced objective flux (FVSEOF) with 
grouping reaction (GR) constraints, in an effort to identify gene amplification targets by considering reactions that 
co-carry flux values based on physiological omics data via "GR constraints". This method scans changes in the 
variabilities of metabolic fluxes in response to an artificially enforced objective flux of product formation. The gene 
amplification targets predicted using this method were validated by comparing the predicted effects with the 
previous experimental results obtained for the production of shikimic acid and putrescine in Escherichia coli. 
Moreover, new gene amplification targets for further enhancing putrescine production were validated through 
experiments involving the overexpression of each identified targeted gene under condition-controlled batch 
cultivation. 

Conclusions: FVSEOF with GR constraints allows identification of gene amplification targets for metabolic 
engineering of microbial strains in order to enhance the production of desired bioproducts. The algorithm was 
validated through the experiments on the enhanced production of putrescine in £ coli, in addition to the 
comparison with the previously reported experimental data. The FVSEOF strategy with GR constraints will be 
generally useful for developing industrially important microbial strains having enhanced capabilities of producing 
chemicals of interest. 

Keywords: Flux variability scanning based on enforced objective flux, Grouping reaction constraints, Putrescine, 
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Background 

One of the most ambitious goals in metabolic engineering 
is the design of biological systems based on in silico pre- 
dictions using mathematical models. The advent of high- 
throughput technologies and the completion of genome 
sequencing for many organisms have led to an explosion 
of systems-wide biological data [1,2]. Genome-scale stoi- 
chiometric models of the increasing number of microor- 
ganisms and mammalian cells have been developed at the 
moment [3,4]. Some of such models have been used to 
identify gene knockout targets for the efficient production 
of important industrial chemicals, including amino acids 
[5,6] and chemicals that are conventionally derived from 
petroleum [7-9]; other such models have been used to 
identify drug targets in pathogens [10-12]. In modeling 
and simulation approaches, target reactions whose knock- 
out is predicted to overproduce the chemical of interest 
can be easily tested experimentally by deleting the corre- 
sponding genes in the microbial host. 

Increasing the expression levels of the relevant genes 
has also been successfully employed for the overproduc- 
tion of target chemicals [13,14]. To avoid unnecessarily 
massive experiments to be performed, several computa- 
tional algorithms have been devised in an effort to reveal 
the relationship between metabolic reactions and the bio- 
logical properties of interest [15-27]; however, the identifi- 
cation of gene amplification targets is more complicated 
than the identification of gene knockout targets; hence, 
correlations among the genes, mRNAs, transcriptional or 
translational regulations, proteins, and metabolic fluxes 
must be carefully examined. Genome-scale metabolic 
models that rely on constraints-based flux analysis with- 
out additional physiological information are limited in 
their ability to describe the complex nature of biological 
systems, particularly biological phenomena beyond me- 
tabolism. Several systematic methods have been devel- 
oped to overcome such limitations: flux variability 
analysis (FVA) [17,19-21], flux coupling analysis [16-18], 
flux sensitivity analysis [15], flux response analysis [26], 
OptReg [22], genetic design through local search [25], 
OptForce [27], and flux scanning based on enforced ob- 
jective flux (FSEOF) [23]. In particular, FSEOF is a 
method that first scans and searches for variations in the 
metabolic fluxes in response to the enforced fluxes direc- 
ted towards a target product. Reactions were then selected 
as amplification targets, the flux values of which increased 
in accordance with the enforced fluxes toward the produc- 
tion of a target chemical. This method was experimentally 
validated by identifying amplification targets that 
improved the production of lycopene in Escherichia coli 
[23]. These approaches demonstrated that incorporating 
physiological constraints during the model simulation are 
critical to identifying trustworthy gene amplification tar- 
gets, but much improvement is still needed [24,28]. One 



of the major problems is the existence of a too large flux 
solution space in optimization problems. 

In this study, in order to systematically handle the large 
flux solution spaces, as also revealed in the implementa- 
tion of FSEOF [23], we considered functionally grouped 
reactions that simultaneously carry fluxes based on unique 
features of microbial genomes. Considering such function- 
ally grouped reactions helps reducing the number of and 
selecting multiple solutions existing for each optimal ob- 
jective value, enabling to identify more reliable gene amp- 
lification targets when combined with FSEOF. Grouped 
reactions were previously revealed by genomic context 
and flux-converging pattern analyses as promising con- 
straints [28]. Genomic context analysis interrogates con- 
served neighborhood, gene fusion, and co-occurrence 
using a STRING database with the goal of suggesting 
groups of reaction fluxes that are most likely correlated in 
their on/off activities [28,29]. Flux-converging pattern ana- 
lysis further limits the range of possible flux values in a 
metabolic reaction by examining the number of carbon 
atoms in metabolites that participate in the reactions and 
the converging patterns of fluxes from a carbon source 
(see Methods and Figure 1) [28]. Consequently, flux bal- 
ance analysis (FBA) with constraints controlling simultan- 
eous on/off activity (C on/q ff) and the flux scale {C sca i e ) of 
the metabolic reactions accurately predicted flux distri- 
butions in gene knockout mutant strains [28]. 

Based on these analyses, the grouping reaction (GR) 
constraints that constrain reactions to co-carry fluxes 
altogether regardless of the condition were incorporated 
into the E. coli genome-scale metabolic model. The model 
then facilitated the scanning of changes in the variability 
among metabolic fluxes using FVA in response to the 
enforced enhancement of the fluxes toward a target chem- 
ical. This newly developed method, called flux variability 
scanning based on enforced objective flux (FVSEOF) with 
GR constraints, was employed in this study to identify 
gene amplification targets for the production of target che- 
micals. FVSEOF with GR constraints was first validated 
based on amplification targets reported for the production 
of shikimic acid and putrescine in E. coli, and then further 
validated by actually engineering E. coli for the enhanced 
production of putrescine based on new amplification 
targets. 

Methods 

E coli genome-scale metabolic model 

EcoMBEL979 was used throughout this study [30], which is 
a slightly modified version of the genome-scale E. coli meta- 
bolic network model, /JR904 [31]. EcoMBEL979 contains 
814 metabolites (144 extracellular metabolites and 670 
intermediates) and 979 metabolic reactions, along with a 
biomass equation derived from the E. coli biomass compo- 
sition [32]. 
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Figure 1 Schematic illustration of the FVSEOF method with GR constraints. Functionally grouped reactions were considered based on 
genomic context and flux-converging pattern analyses obtained from the STRING database. FVSEOF was then performed under GR constraints to 
identify gene amplification candidates for the production of a target chemical. The candidates were evaluated based on the model predictions 
and additional criteria of the flux bias {V avg ) and the slope of the flux changes {q s i ope )- Each rectangle containing a C x J y index and a line with 
different colors defines the reaction groups that are likely on or off simultaneously, as determined by genomic context and flux-converging 
pattern analyses. The C x J y index for each reaction is determined by flux-converging pattern analysis. C x and J y denote the total number of carbon 
atoms in metabolites that participate in each reaction and the type of fluxes through the flux-converging metabolites from a carbon source, 
respectively. The red metabolites indicate flux-converging metabolites. The flux-converging metabolites indicate metabolites at which two 
pathways split by another metabolite recombine. For example, glyceraldehyde-3-phosphate converges the fluxes split by the fructose- 
bisphosphate aldolase from the fructose-6-phosphate. The flux-converging metabolites categorize J y into four types, indicated as J Al J Bl J Q and J D . 
Each subscript of J y denotes the number of flux-converging metabolites that are passed zero, one, two, or three times, respectively, for a given 
flux from a carbon source. The subscript E is specially denoted to indicate the fluxes derived from pyruvate. The values of C x J y for each reaction 
were assigned based on possible flux routes reaching from glucose, and are partitioned by a slash. 



Constraints-based flux analysis 

The stoichiometric relationships among the metabolites and 
the reactions of the E. coli genome-scale metabolic model 
were balanced under the pseudo-steady state assumption. 
The balanced reaction model was almost always underde- 
termined in calculations of the flux distribution due to in- 
sufficient measurements of the extracellular fluxes. Thus, 
the unknown fluxes within the metabolic reaction network 
were calculated by linear programming-based optimization 
using an objective function that maximized the growth rate, 
subject to constraints pertaining to mass conservation and 
reaction thermodynamics [33], This optimization problem 
can be mathematically formulated as follows: 



M 

where Sy represents the stoichiometric coefficient for me- 
tabolite i in reaction v y is the flux of reaction / is the set 
of all reactions, and b t is the net transport flux of metabolite 
L If this metabolite is an intermediate, b t is equal to zero, aj 
and /? y are the lower and upper bounds of the flux of reac- 
tion respectively. Herein, the flux of any irreversible re- 
action is considered to be positive; a negative flux 
indicates the reverse direction of a reaction. 
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Grouping reaction (GR) constraints based on the genomic 
context and flux-converging pattern analyses 

The algorithm introduced in this study, FVSEOF with 
GR constraints, starts with formulation of GR con- 
straints, which are based on the genomic context and 
flux-converging pattern analyses (Figure 1). Briefly, 
genomic context and flux-converging pattern analyses 
aim at grouping functionally related reactions. Such 
functionally related reactions were constrained to be 
on or off simultaneously (Figure 1) [28]. First, reac- 
tions were grouped using STRING database that per- 
forms genomic context analysis, including conserved 
neighborhood, gene fusion, and co-occurrence [28,29]. 
Simultaneous on/off constraint {C on/o ff) can be 
described as follows: 

y(yi) =y(v 2 ) 

y(vi) -ai<vi<y(vi) - p x 

y(v 2 ) • a 2 <v 2 <y(v 2 ) • fi 2 

where y(vj) and y(v 2 ) indicate binary variables (on or 
off) of a certain reaction 1 and 2, respectively. 

Each reaction is then given a CJ y index, determined by 
flux-converging pattern analysis. C x and J y denote the total 
number of carbon atoms in metabolites that participate in 
each reaction and the number of passing flux-converging 
metabolites, respectively. Here, it should be noted that 
cofactors were not considered because the flux scales are 
controlled by the carbon number of primary metabolites, 
not cofactors, according to 13 C-based flux analysis [28]. 
For J yi the flux- converging metabolites indicate metabo- 
lites at which two pathways split by another metabolite 
converge. J y has four types, including J A) J B , J c , and J D , 
depending on the characteristics of flux-converging meta- 
bolites. Subscript of J y denotes the passing number of flux- 
converging metabolites, counting zero, one, two, or three 
times for the flux coming from a carbon source. In some 
cases, the subscript E is placed next to the subscripts of A, 
B, C, or D to indicate the fluxes derived from pyruvate, 
which causes more complex changes in flux distributions. 
The values of CJ y for each reaction were assigned based 
on possible flux routes reaching from glucose, and are 
partitioned by a slash. Based on this analysis, another 
constraint C sca i ei indicating the flux scale of a reaction, 
can be given to the metabolic reactions. First, terms 
used to describe the flux scale of the reaction are as fol- 
lows: 



Cxy/yj 



where C x j indicates the carbon number involved in a re- 
action Jyj the number of the passing of the flux 



through the flux-converging metabolite near reaction 
and Nqr; the total number of carbon of primary meta- 
bolites without cofactors in reaction 

If reaction 1 and 2 were predicted to be in the same 
functional unit according to the genomic context ana- 
lysis, and their C x iJ y \ and C x2 J y2 are equivalent, C scaie is 
applied to these two reactions, which is defined as 
follows: 



-<S 



\ 



where v" and v 2 are the normalized flux of reaction 1 
and 2, obtained by dividing each reaction flux by the 
carbon source uptake rate, such as glucose. 5 is the 
constant defining the flux level of reactions in this 
functional unit; the value of 5 is recommended as 0.3. 

Flux variability scanning based on enforced objective flux 
(FVSEOF) with grouping reaction (GR) constraints 

Once grouping reaction constraints are defined, FVSEOF 
with GR constraints is subsequently performed as follows 
(Figure 2). First, the initial or theoretical minimum 
(C^ product) and theoretical maximum (v^et product) of 
the target product formation rates were calculated; these 
were implemented by setting the objective function as 
minimizing and maximizing the target product formation 
rate using constraints-based flux analysis with GR con- 
straints. This can be formulated as follows: 

Min/Max Z(y target product) = v Targtt product ° r V ?arget product 



Subject to 

k < s ij v j ^ u i 
a j <v j <P j 

Pj = 1000 mmol • g DCW 1 ■ h _1 

a; = -1000 mmol • g DCW" 1 • IT 1 

vro^lOmmol-gDCW-^h- 1 

where ^rget product indicates the initial or minimal point 
of the flux value constrained for the target bioproduct, 
while v^t product indicates the maximal flux value for 
the bioproduct. l t and u t are the lower and upper bound 
for the net transport flux of metabolite z, respectively, 

and Vcarton * s tne car bon source uptake rate. 
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Figure 2 Framework of the FVSEOF with GR constraints for identifying gene amplification targets that enhance the production of a 
target product. The FVSEOF method scans the changes in the metabolic flux variabilities in response to an enhanced flux toward a target 
product. The method then selects amplification target reactions, the fluxes of which increase in response to the forced increase in the flux 
toward the target bioproduct. (A) During the FVSEOF implementation, three types of intracellular flux profiles are typically identified: increased, 
decreased, or unchanged, but oscillatory flux profiles can also be found in some cases. (B) To evaluate the gene amplification candidates, the 
slope (q s io P e) was calculated based on a linear regression between the enforced production rate of a target product and the V avg values of the 
candidate reactions. The positive correlation in the slope indicates that the corresponding reaction may be a gene amplification candidate. On 
the basis of the q sjope values, we considered the sensitivities of the identified gene amplification candidates to the enforced production of a 
target chemical. A large value of q s i ope indicates that the corresponding reaction may be more sensitive to the enforced production of a target 
chemical, than reactions with smaller q sjope values. 



Second, the cell growth rate, Z{v biomas ^ was maxi- 
mized while gradually increasing the target product for- 
mation rate from its initial (or minimal) flux value to its 

near theoretical maximum: v\ 



nforced 



, .initial 



target product target product 



( v max _ ..initial \ jy- Sk\k — 1 9 vi — M 

y target product v target product J > 1X V v l /V x 5^) J 

(n>10) [23]. The v^ d pmduct is an additional constraint 
provided during this stage of the constraints-based flux 
analysis; it starts with the initial value v\ 



initial 



product 



plus one 

n m of the difference between the v?™ etproduct and 



initial 



target product > and is increased to a value adjacent to 
v Z X g et product in /c steps. 

Third, FVA was carried out with GR constraints by 
maximizing or minimizing the fluxes of all intracellular 
reactions, Z(v intraceUular reaction)* with additional con- 
straints: the enforced production rate of the target 



bioproduct, which varied from its initial to maximum 
values in 10 steps, and 95% optimal cell growth rate, 
Vbiomass = 0-95 ■ Z{v biomass ) opt , for each step. The attain- 
able flux ranges of intracellular reactions for each step 
were subsequently subjected to the targeting criteria 
introduced in the following section. 

FVSEOF with GR constraints was calculated using 
mixed integer nonlinear programming with the DICOPT 
solver, subject to the constraints including GR constraints, 
mass conservation and reaction thermodynamics. 

Flux bias, its slope and flux capacity as targeting criteria 

Flux bias (V aV g)> its slope {q s i op e) and flux capacity (l so j) 
were employed as targeting criteria for the initial set of 
gene amplification targets predicted from FVSEOF with 
GR constraints (Figure 2 and 3). Among them, V avg and 
l so i were determined as follows in order to effectively 
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Figure 3 Nine types of changes in the flux patterns based on combinations of positive and negative changes in V avg and / so/ . Types 1, 2, 
and 3 have the amplification candidates positively correlated with the production of a target chemical. Types 4, 5, and 6 were negatively 
correlated with the target product. Types 7, 8, and 9 did not show unique patterns with the enforced fluxes towards the target chemical. Nine 
possible combinations of flux bias {V avg ) and flux capacity (/ so/ ) for each reaction were investigated, and displayed on the bottom of the figure. 



investigate the changes of flux variabilities for genetic 
perturbations [28]: 

V + V . 

-i/ y max 1 y mm 

Vavg ~ 2 

hoi ^max ^min 

The V' max 2Lnd V' m i n indicate the maximal and minimal 
flux values for a reaction under the given condition. The 
I so i indicates the difference between the maximal and 
minimal flux values for a reaction. q s i ope was calculated 
using linear regression of the flux values for a reaction to- 
wards the gradually maximized product formation rate. 

Bacterial strains and plasmids 

The E. coli strains used in this study are listed in the 
Additional file 1. The XQ52 strain, a putrescine produ- 
cer, was used as a base strain [34]. E. coli TOP10 was 



used for gene cloning studies. The plasmid pl5SpeC 
containing a strong tac promoter was used as an expres- 
sion vector. The plasmid pl5SpeC was constructed from 
the pTacl5K plasmid by cloning the speC gene (encod- 
ing ornithine decarboxylase in the putrescine biosyn- 
thetic pathway) into the site between the EcoRI and SacI 
restriction enzyme sites of pTacl5K. The plasmid con- 
tained a kanamycin resistance selective marker. Cells 
were grown in Luria-Bertani (LB) broth or on plates 
containing appropriate antibiotics at 37°C for the con- 
struction of strains and plasmids. Antibiotics were added 
at following concentrations: 50 ug/mL ampicillin, 25 ug/ 
mL kanamycin, and 35 ug/mL chloramphenicol. 

The plasmids used in this study are listed in the Add- 
itional file 1. Polymerase chain reaction (PCR) primers 
for the gene cloning studies conducted here are listed in 
the Additional file 2. Pfu DNA polymerase was pur- 
chased from Solgent (Daejeon, Korea). Restriction 
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enzymes and T4 DNA ligase were obtained from New 
England Biolabs (Ipswich, MA) and Roche (Mannheim, 
Germany), respectively. The genomic DNA of E. coli 
W3110 was amplified to overexpress the target genes 
using the Pfu polymerase and PCR primers (Additional 
file 2). The PCR product was then digested with Sad 
and Xbal, and ligated into pl5SpeC at the same restric- 
tion sites downstream of the tac promoter. 

Fermentation 

Batch cultivation was conducted at 37°C in a 6.6 L jar 
fermentor (Bioflo 3000; New Brunswick Scientific Co., 
Edison, NJ) containing 2 L R/2 medium supplemented 
with 10 g/L glucose and 3 g/L (NH 4 ) 2 S0 4 . The R/2 
medium (pH 6.8) contained (per liter): 2 g (NH 4 ) 2 HP0 4 , 
6.75 g KH 2 P0 4 , 0.85 g citric acid, and 0.7 g 
MgS0 4 -7H 2 0. In addition, 5 mL/L of a trace metal stock 
solution [35] was added. The trace metal solution con- 
tained per liter of 5 M HC1: 10 g FeS0 4 -7H 2 0, 2.25 g 
ZnS0 4 -7H 2 0, 1 g CuS0 4 -5H 2 0, 0.5 g MnS0 4 -5H 2 0, 
0.23 g Na 2 B 4 O 7 10H 2 O, 2 g CaCl 2 -2H 2 0, and 0.1 g 
(NH 4 ) 6 Mo 7 0 24 . One milliliter of the overnight culture 
was transferred into a 300 mL Erlenmeyer flask contain- 
ing 50 mL of the R/2 medium at 37°C and spun at 220 
rpm in a shaking incubator (JEIOTech. Co. SI-900R). 
After obtaining an initial OD 600 of 0.3, the seed cultures 
(200 mL) were introduced into the bioreactor for batch 
cultivation. The culture pH was maintained at 6.8 by the 
addition of 6 M KOH. The dissolved oxygen concentra- 
tion was maintained at 20% air saturation by automatic- 
ally adjusting the agitation speed. Under the comparable 
batch culture conditions, the single gene-overexpressing 
strains based on the E. coli XQ52 strain harboring 
pl5SpeC, denoted as XQ52 (pl5SpeC), with each target 
gene were tested by flask cultivation in duplicate using 
R/2 medium supplemented with 10 g/L glucose at 37 °C. 

Analytical procedures 

Cell growth was estimated by measuring the optical density 
at 600 nm (OD 600 ) using an Ultrospec 3000 spectropho- 
tometer (Amersham Biosciences, Uppsala, Sweden). Glu- 
cose concentrations were measured using a glucose 
analyzer (model 2700 STAT; Yellow Springs Instrument, 
Yellow Springs, OH, USA). The concentrations of glucose 
and organic acids were determined by high-performance li- 
quid chromatography (ProStar 210; Varian, Palo Alto, CA) 
equipped with UV/visible light (ProStar 320; Varian, Palo 
Alto, CA) and refractive index (Shodex RI-71, Tokyo, 
Japan) detectors. A MetaCarb 87H column (300 by 7.8 
mm; Varian) was eluted isocratically with 0.01 NH 2 S0 4 at 
60°C at a flow rate of 0.4 mL/min. 

The putrescine concentration was determined by deri- 
vatizing putrescine with o-phthaldialdehyde (OPA; 
Sigma, St. Louis, MO), and the o-phthaldialdehyde 



derivative was detected by high-performance liquid 
chromatography (1100 Series HPLC, Agilent Technolo- 
gies, Palo Alto, CA) with UV detection, as described pre- 
viously [34]. The OPA derivatization reagent was 
prepared as described previously [34,36,37]. Following 
the addition of the OPA reagent, the mixture was fil- 
tered through a 0.2 mm PVDF syringe filter (Whatman, 
Maidstone, UK), and the filtrate was immediately 
injected into the HPLC. A SUPELCO C18 column (cat# 
504955; 5|im, 150 mm x 4.6 mm) was operated at 25°C 
with a 0.8 mL/min mobile phase flow rate. The mobile 
phase consisted of solution A (55% methanol in 0.1 M 
sodium acetate, pH 7.2) and solution B (methanol). The 
following gradient was applied (values given in vol%): 1- 
6 min, 100% A; 6-10 min, linear gradient of B from 0% 
to 30%; 10-15 min, linear gradient of B from 30% to 
50%; 15-19 min, linear gradient of B from 50% to 100%; 
19-23 min, 100% B; 23-25 min, linear gradient of B 
from 100% to 30%; 25-28 min, linear gradient of B from 
30% to 0% [34]. The derivatized putrescine was detected 
at a wavelength of 230 nm using a variable wavelength 
detector (G1314A, Agilent Technologies). 

Results and discussion 

FVSEOF with GR constraints 

Functionally related reactions can be grouped by gen- 
omic context and flux-converging pattern analyses [28]. 
Several reactions appeared to be related with one an- 
other based on genomic context analysis of conserved 
neighborhoods, gene fusions, and co-occurrence [28,29]. 
Flux-converging pattern analysis narrows the range of 
plausible flux values for metabolic reactions by examin- 
ing the number of carbon atoms in metabolites that par- 
ticipate in reactions and the converging patterns of 
fluxes from a carbon source [28]. By controlling the sim- 
ultaneous on/off activity {C on/o ff) and flux scale (C sca i e ) 
of the metabolic reactions, based on FBA with GR con- 
straints, the flux distributions in gene knockout mutants 
were accurately predicted [28]. In this study, the GR 
constraints were further applied to reactions related to 
the biosynthesis of a target chemical to improve the 
model accuracy (Figure 1 and the Additional file 3). 

FVSEOF with GR constraints was implemented as fol- 
lows (Figure 2). First, the theoretical minimal and max- 
imal flux values for the target product formation were 
calculated using constraints -based flux analysis by min- 
imizing and maximizing the target product formation 
rate with GR constraints. Second, again with GR con- 
straints, the cell growth rate was maximized while grad- 
ually increasing the constraint value for the target 
product formation rate (our objective function which is 
artificially enforced) from a minimal to a theoretical 
maximum, as calculated from the first step. Finally, FVA 
was conducted by maximizing and minimizing the fluxes 
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of all intracellular reactions under additional constraints, 
including GR constraints, the enforced production rate 
of the target chemical varied from a minimal to a max- 
imal value, and a 95% optimal growth rate constraint for 
each step. The attainable flux ranges for the intracellular 
reactions were calculated under the imposed constraints 
for each of the three steps. 

Criteria for selecting gene amplification targets 

Initial simulation results for FVSEOF with GR con- 
straints were filtered based on rational criteria in an ef- 
fort to select only the most effective amplification 
targets. The most important criterion was to identify 
gene amplification targets, the fluxes of which increased 
with the flux directed toward the target chemical. This 
procedure was implemented with quantitative values of 
q s io P e> Vavg, and I (Figures 2 and 3). The flux bias (V avg ) 
and flux capacity (l so i) indicate an average value for the 
maximal and minimal flux values and the length of the 
attainable flux ranges for a reaction, respectively [28]. Fi- 
nally, the gene amplification candidates were evaluated 
by calculating the slope (q s iope) of the V avg flux for each 
metabolic reaction using linear regression analysis 
(Figure 2). Changes in the patterns of the reaction fluxes 
in response to incrementally increasing fluxes toward a 
target product were categorized into nine types based on 
combinations of positive and negative changes in V avg 
and l so i in order to facilitate the identification of amplifi- 
cation targets (Figure 3). Types 1, 2, and 3 displayed 
positive correlations with the amplification candidates 
for the production of a target chemical; types 4, 5, and 
6 displayed negative correlations with the amplification 
candidates. Finally, types 7, 8, and 9 displayed no clear 
correlations with the production of a target chemical, 
based on V avg (Figure 3). The reaction sets that were 
positively correlated (type 1, 2, and 3) were initially 
selected as amplification candidates. 

Reaction sets belonging to type 1, 2 and 3, showing 
positive correlations with the enforced fluxes toward a 
target chemical, can then be further divided into strongly 
and weakly positive reactions (Figure 3 and Additional 
file 4). This step also allows narrowing down the candi- 
dates of gene amplification targets. The strongly positive 
reactions display a continuously increasing V avg and a 
positive q s io P e in response to the enhanced production of 
a target chemical, whereas weakly positive reactions 
show the same pattern, except for the presence of a par- 
tially negative q s i op e (Additional file 4). Certainly, 
strongly positive reactions deserve primary attention as 
potential gene amplification targets. 

The potential gene amplification candidates were 
prioritized by considering the l so i value, which indicates 
the length between the maximal and minimal flux values 
of a metabolic reaction. Among the reactions that were 



positively correlated with the desired product, reactions 
with smaller values of l so i received higher priorities be- 
cause these reactions were more likely to display the 
predicted flux values than reactions with larger values of 
l so i. A final list of gene amplification targets obtained 
from the above procedure was then selected based on 
biological knowledge. 

Implementation of FVSEOF with GR constraints for 
shikimic acid production in E. coli 

FVSEOF with GR constraints was employed to identify 
gene amplification targets for the enhanced production 
of an important aromatic chemical shikimic acid in E. 
coli (Figure 2). Shikimic acid is a key metabolic inter- 
mediate in the aromatic amino acid biosynthetic path- 
way. Shikimic acid and its derivatives are industrially 
important starting compounds for the production of sev- 
eral chemicals, such as phenols, herbicides, antibacterial 
agents, and the neuramidase inhibitor Tamiflu used for 
the treatment of influenza infections [38,39]. FVSEOF 
with GR constraints predicted that 11 reaction fluxes in 
the glycolysis (glk and pps), pentose phosphate pathway 
(rpi, talAB, and tktAB), and the shikimic acid biosyn- 
thetic pathway (aroB, aroD, aroE, aroF, aroG, and aroH) 
were potential amplification targets. The amplification of 
aroB, aroD, aroE, aroF, aroG, aroH, talAB, tktA, glk, and 
pps genes [38-46], which are the amplification targets 
predicted by FVSEOF with GR constraints, was previ- 
ously reported to enhance the production of shikimic 
acid. The previous FSEOF method without FVA and GR 
constraints could not identify pps gene as an amplifica- 
tion target. FSEOF results without FVA and GR con- 
straints did not show notable fluxes among metabolic 
reactions controlled by the pps gene in response to the 
enforced shikimic acid production rate; however, the 
FVSEOF method with GR constraints correctly pre- 
dicted the pps gene as one of the amplification targets 
beneficial for the accumulation of phosphoenolpyruvate, 
an important precursor for the production of shikimic 
acid from pyruvate. This consistency partly demon- 
strated the power of utilizing FVA and GR constraints 
for predicting reliable amplification targets by FSEOF. In 
practice, the overexpression of phosphoenolpyruvate 
synthase encoded by the pps gene also increased the 
yield of precursors for the production of shikimic acid 
[43,45] . Thus, this strategy enabled the successful identi- 
fication of gene amplification targets for the enhanced 
production of a primary metabolite, shikimic acid, in E. 
coli, in accordance with previous literature reports. 

Implementation of FVSEOF with GR constraints for 
enhanced putrescine production in E coli 

The general applicability of FVSEOF with GR con- 
straints was examined by applying the method to 
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putrescine production in E. coli. Putrescine (1,4-diami- 
nobutane) is an important industrial precursor for the 
synthesis of polymers, pharmaceuticals, surfactants, and 
certain additives [34]. We confirmed the validity of the 
newly predicted gene amplification targets by compari- 
son with the genes engineered in the previously reported 
putrescine-producing E. coli XQ52 (pl5SpeC) strain 
[34]. FVSEOF with GR constraints predicted potential 
gene amplification targets among the reactions involved 
in glycolysis (eno, pgm, gapA, fbaAB, tpiA, pgk, pyI<AF, 
and glk), TCA cycle (icd, acnA, acnB, and gltA), putres- 
cine biosynthesis (gdhA, argA, argB, argC, argD, argE, 
speC, and speF), and other pathways [ackA and ppc). 
The predicted amplification targets {argB, argC, argD, 
argE, speC, and speF) involved in the putrescine biosyn- 
thetic pathway were consistent with the mutations intro- 
duced in the E. coli XQ52 (pl5SpeC) strain, as described 
in the previous report [34]. 

The genes predicted to be relevant to the putrescine 
biosynthetic pathway were expected based on the 



pathway knowledge, and were intuitively obvious; hence, 
we focused on the effects associated with amplifying the 
predicted gene targets involved in other metabolic path- 
ways in order to more rigorously validate FVSEOF with 
GR constraints. Accordingly, each of the predicted amp- 
lification targets was examined one by one by amplifying 
the gene dosage in the E. coli XQ52 (pl5SpeC) strain 
(see Methods). Among these genes, the glk, acnA, acnB, 
ackA, and ppc gene, five out of the sixteen target genes, 
were found to attain increased putrescine yield when 
they were individually amplified in the E. coli XQ52 
(pl5SpeC) strain (Additional file 5). These strains, ini- 
tially examined using flask cultivation, were further vali- 
dated by batch cultivation at 37°C under aerobic 
condition (Figure 4 and see the Additional files 1 and 2). 
The recombinant E. coli XQ52 (pl5SpeC) strain add- 
itionally expressing the glk, acnA, acnB, ackA, and ppc 
genes resulted in the production of 2.23, 1.90, 1.89, 2.04, 
and 2.06 g/L of putrescine, respectively, which are 20.5 
% more than that (1.68 g/L) produced by the control 
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Figure 4 Batch cultivation of the control and engineered E. coli strains producing putrescine. (A) XQ52 (p15SpeC) (control), (B) XQ52 
(p15SpeC-Glk), (C) XQ52 (p15SpeC-AcnA), (D) XQ52 (p15SpeC-AcnB), (E) XQ52 (p15SpeC-AckA), and (F) XQ52 (p15SpeC-Ppc) strains were tested 
for putrescine production. The symbols indicate: cell concentration measured by OD 600 (•), glucose concentration (■), and putrescine 
concentration (A) for the engineered strains in (A) to (F), and putrescine concentration (A) for the XQ52 (p15SpeC) (control) strain in (A). 



Park et a I. BMC Systems Biology 201 2, 6:1 06 
http://www.biomedcentral.eom/1 752-0509/6/1 06 



Page 10 of 11 



strain on average (Figure 4). The yields of putrescine 
obtained with these strains were 0.223, 0.190, 0.189, 
0.204, and 0.206 g putrescine per g glucose, which are 
again higher than that (0.168 g putrescine per g glucose) 
obtained with the control strain. Thus, FVSEOF with 
GR constraints could be successfully used to identify 
non-obvious gene amplification targets that enhance the 
production of putrescine in E. coli. 

Other gene amplification targets identified by FVSEOF 
with GR constraints, which did not affect putrescine 
production in flask cultivation, also deserve discussion. 
These false-positive hits are most likely involved in bio- 
logical processes that were not accurately captured in 
the genome-scale metabolic model. These ineffective 
genes, including eno, pgm, gap A, fbaAB, tpiA, pgk, and 
pykAE genes in glycolysis, and icd and git A genes in 
TCA cycle might have been associated with transcrip- 
tional and translational regulations because the direct 
correlation between gene expressions and the metabolic 
fluxes was not observed. The fact that some of obvious 
gene amplification targets, such as icd gene responsible 
for biosynthesis of a-ketoglutarate, seem to be resistant 
to gene manipulations indicates that other biological 
variables may affect the effects of the gene amplifica- 
tions. Potential variables include the plasmid copy num- 
ber, gene dosage, optimal gene expression, and the gene 
expression method, either plasmid-based overexpression 
or chromosomal integration [47]. Although we improved 
the accuracy of the predicted gene targets by imposing 
GR constraints, these factors should be carefully consid- 
ered in any implementation of the FVSEOF method with 
GR constraints [48]. 

Conclusions 

FVSEOF with GR constraints, which is an upgraded ver- 
sion of the FSEOF method, allows for the in silico identi- 
fication of fluxes to be amplified for the enhanced 
production of target products. This method was con- 
ducted through the analysis of trends in reaction flux 
variability in response to varying the flux of target chem- 
ical production from initial to maximal flux values under 
GR constraints. The confidence with which amplification 
targets are identified may be increased by incorporating 
physiological data. This approach involves grouping 
functionally related reactions based on their genomic 
context and flux- converging pattern analyses. The inter- 
action data may be obtained easily from public data- 
bases, and subjected to GR constraints. FVA was also 
performed to overcome the problems associated with 
multiple solutions for an optimal objective value. 
FVSEOF with GR constraints was shown to suggest suc- 
cessful metabolic engineering strategies (in particular, 
gene amplification) for the production of shikimic acid 
and putrescine in E. coli. In conclusion, the strategy 



reported here should be generally useful for developing 
industrial strains that display enhanced production of a 
target chemical. 
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